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^^ ■ Black holes are popping up all over the place: in compact binary X-ray sources and GRBs, 

^ ' in quasars, AGNs and the cores of all bulge galaxies, in binary black holes and binary black 

O ■ hole-neutron stars, and maybe even in the LHC! Black holes are strong-field objects governed 

1/'^ ' by Einstein's equations of general relativity. Hence general relativistic, numerical simulations 

of dynamical phenomena involving black holes may help reveal ways in which black holes can 
form, grow and be detected in the universe. To convey the state-of-the art, we summarize several 
representative simulations here, including the collapse of a hypermassive neutron star to a black 
ri [ hole following the merger of a binary neutron star, the magnetorotational collapse of a massive 

^ |. star to a black hole, and the formation and growth of supermassive black hole seeds by relativistic 

' ' MHD accretion in the early universe. 

1. Introduction 

Black holes are 'sighted' everywhere in the universe these days. Originally located 
^ ' in compact binary X-ray sources in the 1970's, the cosmic presence of black holes has 

f — I expanded considerably in recent decades. They now are believed to be the central engines 

CO . that power quasars, active galactic nuclei (AGNs) and gamma-ray bursts (GRBs). They 

are identified in the cores of all bulge galaxies. They are presumed to form signficant 
populations of compact binaries, including black hole-black hole binaries (BHBHs) and 
black hole- neutron star binaries (BHNSs). Black holes may even show up soon in the 

l^^_^ ■ Large Hadron Collider! 

f^ I Gravitationally, black holes are strong-field objects whose properties are governed by 

Einstein's theory of relativistic gravitation — general relativity. General relativistic sim- 
ulations of gravitational collapse to black holes, BHBH mergers and recoil, black hole 

S^ ■ accretion, and other astrophysical phenomena involving black holes may help reveal how, 

H I when and where black holes form, grow and interact in the physical universe. As a con- 

sequence, such simulations can help identify the ways in which black holes can best be 
detected. 

To illustrate how our understanding of black hole phenomena is sharpend by large- 
scale simulations in general relativity, we summarize in this paper the results of several 
recent computational investigations. The first few involve the formation of black holes 
from stellar collapse, while the last one concerns supermassive black hole growth via 
disk accretion. Most of these simulations utilize the tools of numerical relativity to solve 
Einstein's field equations of general relativity. So we shall begin with a brief overview of 
this important, rapidly maturing field. 

2. Numerical Relativity and the 3 + 1 Formalism 

Most current work in numerical relativity is performed within the framework of the 
3-1-1 decomposition of Einstein's field equations using some adaptation of the standard 
ADM equations ( [Arnowitt, Deser fc Misner 1962D . In this framework spacetime is sliced 
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Figure 1. 3+1 decomposition of spacetime. 



up into a sequence of spacelike hypersurfaces of constant time t, appropriate for solving 
an initial-value problem. Consider two such time slices separated by an infinitesimal 
interval dt as shown in Fig [1] The spacetime metric measures the invariant interval 
between neighboring points A and B on the two slices according to 



ds^ 



-a^dt^ + y,j{dx' + P'dt){dx^ + (3^ dt) 



(2.1) 



where 7jj is the spatial 3-metric on a time slice, a is the lapse function determining the 
proper time between the slices as measured by a time-like normal observer n" at rest in 
the slice, and /?' is the shift vector, a spatial 3- vector that describes the relabeling of the 
spatial coordinates of points in the slice. The gravitational field satisfies the Hamiltonian 
and momentum constraint equations on each time slice, including the initial slice at t = 0: 



R + K^ ~ K,jK'^ = 167rp (Hamihonian) 

Dj{K'' - -f'^K) == 8t:S' (momentum) . 



(2.2) 
(2.3) 



Here R = R\ is the scalar curvature on the slice, Rij is the 3-Ricci tensor, Kij is the 
extrinsic curvature, K is its trace, and Dj is the covariant derivative operator on the slice. 
The quantities p and S"* are the mass and momentum densities of the matter, respectively; 
such matter source terms are formed by taking suitable projections of the matter stress- 
energy tensor T"^'' with respect to the normal observer. Included in this stress-energy 
tensor are contributions from all the nongravitational sources of mass-energy, which we 
are simply calling the "matter" (e.g., baryons, electromagnetic fields, neutrinos, etc.). 

A gravitational field satisfying the constraint equations on the initial slice can be 
determined at future times by integrating the evolution equations, 



daij = -2aK,j + D,f3j + Dj(3i 



dtKij = a{Rij 



P'^dkK^, 



2K^kK''^ + KK,j 



K.kdjP'' 



DiDja 



(2.4) 
(2.5) 



-s72]{S - P)) 



where S and 5"^^ are additional matter source terms. The evolution equations guarantee 
that the field equations will automatically satisfy the constraints on all future time slices 
identically, provided they satisfy them on the initial slice. Of course, this statement ap- 
plies to the analytic set of equations and not necessarily to their numerical counterparts. 
Note that the 3-1-1 formalism prescribes no equations for a and (3^ . These four functions 
embody the four-fold gauge (coordinate) freedom inherent in general relativity. Choosing 
them judiciously, especially in the presence of black holes, is one of the main challenges 
of numerical relativity. 
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2.1. The BSSN Scheme 

During the past decade, significant improvement in our abihty to numerically integrate 
Einstein's equations stably in full 3+1 dimensions has been achieved by recasting the 
original ADM system of equations. One such reformulation is the so-called BSSN scheme 
(jShibata fc Nakamura 19951 [Baumgarte fc: Shapiro 1999[ ). In this scheme, the physical 
metric and extrinsic curvature variables are replaced in favor of the conformal metric and 
extrinsic curvature, in the spirit of the "York-Lichnerowicz" split (jLichnerowicz 19441 
lYork 19711) : 

7^j-=e-4Sj-, where e^"^ = 71/^ , (2.6) 

A^3 = k,, - ^%K . (2.7) 

Here a tilde " denotes a conformal quantity and 7 is the determinant of 7,,, . At the same 
time, a connection function F* is introduced according to 

f' = 7^'=r,fe = ~d,f' . (2.8) 

The quantities that are independently evolved in this scheme are now jij , Aij , ^, K and 
F'. The advantage is that the Riemann operator appearing in the evolution equations 
(cf. eqn. (|2.5p ) takes on the form, 

R^J=-l¥"'^^^l%+ 7fe(,aj)P +••• . (2.9) 

' Laplacian' remaining 2nd dcrivs 

Thus the principal part of this operator, ^''"'' dmdijij is that of a Laplacian acting on the 
components of the metric jij . All the other second derivatives of the metric have been 
absorbed in the derivatives of the connection functions. The coupled evolution equations 
for 7y and Aij (cf. eqns. (|2.4p and (|2.5p then reduce essentially to a wave equation, 

d^j,, ~ dtAj ^ i?y ^ V27y . (2.10) 

Wave equations not only reflect the hyperbolic nature of general relativity, but can be 
implemented numerically in a straight-forward and stable manner. By now, numerous 
simulations have demonstrated the dramatically improved stability achieved in the BSSN 
scheme over the standard ADM equations, and considerable effort has gone into explain- 



ing the improvement on theoretical grounds [see, e.g., Baumgarte fc Shapiro 2003 for 



discussion and references] . Many of the recent BHBH merger calculations have been per- 



formed using this scheme, beginning with Campanelli et al. (2006) and Baker et al. (2006) 



(but see lPretorius 2005l for an alternative approach). The same is true for the simulations 
described below. 



3. Binary Neutron Star Mergers and Hypermassive Stars 

The protagonist of several different astrophysical scenarios probed by recent numerical 
simulations is a hypermassive star, typically a hypermassive neutron star (HMNS). A 
hypermassive star is an equilibrium fluid configuration that supports itself against grav- 
itational collapse by differential rotation. Uniform rotation can increase the maximum 
mass of a nonrotating, spherical equilibrium star by at most ^ 20%, but differential rota- 
tion can achieve a much hi gher increase ([Baumgarte, Shapiro fc Shibat a 2000[ Morrison, 
Baumgarte & Shapiro 2004). Dynamical simulations using the BSSN scheme demonstrate 
( [Baumgarte, Shapiro fc Shibata 2000] ) that hypermassive stars can be constructed that 
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are dynamically stable, provided the ratio of rotational kinetic to gravitational potential 
energy, /3, is not too large; for (3 > 0.24 the configuration is subject to a nonaxisymmetric 
dynamical bar instability ( [Shibata, Baumgarte fc Shapiro 2000{|Saijo et al. 200ip . How- 
ever, all hypermassive stars are secularly unstable to the redistribution of angular angular 
momentum by viscosity, magnetic braking, gravitational radiation, or any other agent 
that dissipates internal shear. Such a redistribution tends to drive a hypermassive star to 
uniform rotation, which cannot support the mass against collapse. Hence hypermassive 
stars are transient phenomena. Their formation following, for example, a NSNS merger, 
or core collapse in a massive, rotating star, may ultimately lead to a 'delayed' collapse 
to a black hole on secular (dissipative) timescales. Such a collapse will be accompanied 
inevitably by a delayed gravitional wave burst ( Baumgarte, Shapiro fc Shibata 2000] 
[Shapiro 2000| ). 

The above scenario has become very relevant in light of the most recent and detailed 
simulations in full general relativity of NSNS mergers. State-of-the-art, fully relativistic 
simulations of NSNSs have been performed by Shibata and his collaborators (Shibata, 
Taniguchi & Uryu 2003.2005: [Shibata 20051 [Shibata fc Taniguchi 20"06l). They consider 



mergers of n = 1 polytropes, as well as configurations obeying a more realistic nuclear 
equation of state (EOS). They treat mass ratios Qm hi the range 0.9 ^ Qm ^ Ij 
consistent with the range of Qm in observed binary pulsars with accurately determined 



masses ( Thorsett fc Chakrabarty 1999j Stairs 2004)) . The key result is that there exists 



a critical mass Merit ~ 2.5 — 2.7Mq of the binary system above which the merger leads 
to prompt collapse to a black hole, and below which the merger forms a hypermassive 
remnant. With the adopted EOS, the HMNS remnant undergoes delayed collapse in 
about ^ 100 ms and emits a delayed gravitational wave burst. Most interesting, prior to 
collapse, the remnant forms a triaxial bar when a realistic EOS is adopted (see Fig. d]) 
and the bar emits quasiperiodic gravitational waves at a frequency / '^ 3 — 4 kHz. Such a 
signal may be detectable by Advanced LIGO. It is interesting that for the adopted EOS, 
the mass Merit is close to the value of the total mass found in each of the observed binary 
pulsars. Given that the masses of the individual stars in a binary can be determined 
by measuring the gravitational wave signal emitted during the adiabatic, inspiral epoch 
prior to plunge and merger, the detection (or absence) of any quasiperiodic emission from 
the hypermassive remnant prior to delayed collapse may significantly constrain models 
of the nuclear EOS. 

The possibility that a HMNS remnant forms following a NSNS merger had been fore- 
shadowed in earli er Newtonian simulations (Rasio fc Shapiro 1994,1999; Zhuge, Centrella 
& McMillan 1996), in post-Newtonian simulations (Faber fc Rasio 2000,2002) and in con- 
formally flat general relativistic simulations ( [Faber, Grandclement fc Rasio 2004] ). How- 
ever, the recent fully relativistic simulations by Shibata, Taniguchi fc Uryu (2003,2005), 
Shibata (2005) and Shibata fc Taniguchi (2006) provide the strongest theoretical evidence 
of this phenomenon to date, although the details undoubtedly depend on the adopted 
EOS. Triaxial equilibria can arise only in stars that can support sufficiently high values 
of (3 exceeding the classical bifurcation point at /3 « 0.14; reaching such high values 
requires EOSs with adiabatic indicies exceeding F ^ 2.25 in Newtonian configurations, 
and comparable values in relativistic stars. It is not yet known whether the true nuclear 
EOS in neutron stars is this stiff, or what agent for redistributing angular momentum in 
a hypermassive star dominates (e.g., magnetic fields, turbulent viscosity or gravitational 
radiation) . But it is already evident that hypermassive stars are likely to form from some 
mergers and that they will survive many dynamical timescales before undergoing de- 
layed collapse. The recent measurement ([Nice et al. 2005^ of the mass of a neutron star 
in a neutron star-white dwarf binary of M = 2.1Mq establishes an observational lower 
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Figure 2. Formation of triaxial HMNS remnant following NSNS merger in 2.7 Mq system. 
Snapshots of density contours are shown in the equatorial plane. The number in the upper 
left-hand corner denotes the elapsed time in ms; the initial orbital period is 2.11 ms. Vectors 
indicate the local velocity field. [From Shibata, Taniguchi and Uryti 2005.] 



limit for the maximum mass of a neutron star; such a high value suggests that merg- 
ers in typical NSNSs may form hypermassive stars more often than undergoing prompt 
collapse. 



3.0.1. Collapse of a Magnetized HMNS 

The HMNSs found above may survive for many rotation periods. However, as we 
stated, on longer timescales magnetic fields will transport angular momentum and this 
will trigger gravitational collapse. Two important magnetohydrodynamic (MHD) mech- 
anisms which tr ansport angular momentum are magnetic braking ( Baumgarte, Shapiro 
& Shibata 2000, [Shapiro 2000[ |Cook, S hapiro fc Stephens 2003', 'Liu fc Shapi ro 20031 ) 
and the magnetorotational instability (MRI; Vclikhov 1959, Chandrasckhar 196Q1 Balbus 
& Hawley 1991,1998). Magnetic breaking transports angular momentum on the Alfven 
time scale, ta ~ R/va ~ 1(^/10^'^ G)'^ s, where R is the radius of the HMNS. MRI 
occurs wherever angular velocity f2 decreases with cylindrical radius w. This instability 
grows exponentially with an e- folding time of tmri ~ 4 (951/9 In tn)" , independent of 
the field strength. For typical HMNSs considered here, tmri ^ 1 ms. The length scale 
of the fastest growing unstable MRI modes, Amri, does depend on the field strength: 
Amri ~ 3 m (rj/4000s-i)-i {B/W^G) < R. When the MRI saturates, turbulence con- 
sisting of small-scale eddies often develops, leading to angular momentum transport on 
a timescale much longer than tmri ■ The computational challenge of evolving an HMNS 
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is having sufficient spatial grid to resolve the MRI wavelength and sufficient integration 
time to follow the evolution on the long Alfven timescale. 

To determine the final fate of the HMNS, it is necessary to carry out MHD simulations 
in full general relativity. Such simulations have only recently become possible. Duez et al. 
(2005) and Shibata & Sekiguchi (2005) have developed new codes to evolve the coupled 
set of Einstein-Maxwell-MHD equations self-consistently. Our two codes have since been 
used to simulate the evolution of magnetized HMNSs (Duez et al. 2006a, 2006b), and 
implications for short GRBs have been investigated (jShibata et al. 2006) . Both codes 
give very similar results. 

We assume axisymmetry and equatorial symmetry in all of these simulations. We 
use uniform computational grids with sizes up to 500x500 spatial zones in cylindrical 
coordinates. To model the remnant formed in binary merger simulations, we use as our 
initial data an equilibrium HMNS constructed from a F = 2 polytrope with mass M 1.7 
times the spherical mass limit, or 1.5 times the limit for a uniformly rotating star built 
out of the same EOS. These are the typical values expected for a HMNS formed from a 
NSNS merger, for which M r^ 2 x I.AMq = 2.8Mq. The differential rotation proffie is 
chosen so that the ratio of equatorial to central fi is ~ 1/3, comparable to what is found 
in simulations of NSNS mergers. (We find that an HMNS with a more realistic hybrid 
EOS rather than a polytrope evolves similarlv: IDuez et al. 2006bllShibata et al. 2006[) We 
add a seed poloidal magnetic field with strength proportional to the gas pressure. The 
initial magnetic pressure is set much smaller than the gas pressure, but not so small that 
Amri cannot be resolved. Therefore, we set Amri ~ i?/10, corresponding to i? « 10^^ 
G and max(i?^/P) ~ 10~^. The resulting Alfven timescale is about 16 central rotation 
periods, or 600M. 

In our evolutions, the effects of magnetic winding are reflected in the generation of a 
toroidal B field which grows linearly with time during the early phase of the evolution, 
and saturates on the Alfven timescale. The effects of MRI are observed in an exponential 
growth of the poloidal field on the Amri scale, a growth which saturates after a few rota- 
tion periods. The magnetic fields cause angular momentum to be transported outward, so 
that the core of the star contracts while the outer layers expand. After about 66 rotation 
periods, the core collapses to a black hole. Using the technique of black hole excision (see 



Duez, Shapiro & Yo 2004 and references therein) to remove the interior of the black hole, 
with its nasty spacetime singularity, and replace it with suitable boundary conditions on 
all the variables just inside the horizon, we continue the evolution to a quasi-stationary 
state. The final state consists of a black hole of irreducible mass ^ 0.9M surrounded by 
a hot accreting torus with rest mass ^ O.IM and a magnetic field coUimated along the 
rotation axis; see Fig. [3l At its final accretion rate, the torus should survive for ~10ms. 
The torus is optically thick to neutrinos, and we estimate that it will emit ~ lO^^ergs in 
neutrinos before being accreted. We also find that the cone outside the black hole cen- 
tered along the rotation axis is very baryon-poor. All these properties make this system 
a promising central engine for a short-hard GRB. 



4. Black Holes as Central Engines for GRBs 

The combined observations of BATSE, Swift, HETE-2, Chandra and the HST indi- 
cate that GRBs comprise at least two classes: long-soft and short-hard. Long-soft GRBs 
have characteristic timescales r in the range r ~ 2 — lO'^s. They are found in star- 
forming regions (spirals) and some are observed to be associated with supernovae. The 
favored model for the progenitor of a long-soft GRB is the collapse of a massive, rotating, 
magnetized star to a black hole ('collapsar'; |MacFadyen fc: Woosley 1999[ ). By contrast. 
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Figure 3. Collapse of a magnetized HMNS to a black hole. The upper 3 panels show 
snapshots of the rest-mass density contours and velocity vectors on the meridional plane. 
The lower panels show the field lines (lines of constant vector potential A^) for the 
poloidal magnetic field at the same times as the upper panels. The density contours 
are drawn for p/pmax,o = lO"'''^'"'''''^ (j = 0-12). The field lines are drawn for 
A^ = A0,min + (^^,max — A0,min)i/2O (i = 1-19), whcrc A,/,, max aud A^^miu are the maxi- 
mum and minimum value of A^ respectively at the given time. The thick solid curves in the 
lower left corner denote the apparent horizon. [From lDuez et al. 2006al ] 

short-hard GRBs have characteristic timescales in the range r ~ 10ms —2s. They are 
identified in low star-forming regions (ellipticals) where associations with supernovae can 
be excluded. The favored model for their progenitors are either NSNS or BHNS mergers. 
Alternative routes by which NSNS mergers can result in the generation of a short-hard 
GRB are traced in Fig. [H These alternatives have emerged from detailed simulations in 
general relativity. The HMNS route has already been summarized in Section [3l 

The inspiral and merger of NSNSs and BHNSs have important implications for the de- 
tection of gravitational waves with Advanced LIGO. Recent estimates for the rates for de- 
tectable N SNS mergers are promising, in the neighborhood of 20— 30yr~^ ( O'Shaughnessy 
et al. 2006). Simulations in general relativity now underway should be helpful in prepar- 
ing for the exciting possibility of the simultaneous detection of a gravitational wave and 
GRB from the same source in the near future. 



5. Magnetorotational Collapse of Massive Stars to Black Holes 

Recently, we performed simulations in axisymmetry of the magnetorotational collapse 
of very massive stars in full general relativity ( |Liu, Shapiro fc Stephens 2007[ ). Our 
simulations are directly applicable to the collapse of supermassive stars with masses 
M > 10'^ Mq and to very massive Population III stars. They are also relevant for core 
collapse in massive Population I stars, since in all of these cases the governing EOS 
up to the appearance of a black hole can be approximated by an adiabatic T = 4/3 
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Figure 4. Plausible routes for the formation of a short-hard gamma-ray burst (SGRB) central 
engine following the merger of a binary neutron star. Here "GW" ( "B-field J-transport" ) denotes 
angulaj momentum dissipation dominated by gravitational wave emission (magnetic fields). 
[From IShibata fc Taniguchi 2006] ] 



law (although its physical origin is different). These simulations may help explain the 
formation of th e central engine in the coUapsar model of long-soft GRBs ( MacFadyen 
& Woosley 1999). Moreover, some long-soft GRBs observed at very high redshift might 
be related to t he gravitational collapse of very massive Pop III stars ( Schneider, Guetta 
& Ferrara 2002. IBromm fc Loeb 2006^ . Hence these simulations may also provide direct 
insights into the formation of GRB central engines arising from first-generation stars. 

The simulations of Liu, Shapiro fc Stephens (2007) model the initial configurations 
by 71 = 3 polytropes, uniformly rotating near the mass-shedding limit and at the onset 
of radial instability to collapse. These simulations extend the earlier results of Shibata 
& Shapiro 2002 by incorporating the effects of a magnetic field and by tracking the 
evolution for a much longer time after the appearance of a central black hole. The ratio 
of magnetic to rotational kinetic energy in the initial stars is chosen to be small (1% and 
10%). We find that such magnetic fields do not affect the initial collapse significantly. 
The core collapses to a black hole, after which black hole excision is employed to continue 
the evolution long enough for the hole reach a quasi-stationary state. We find that the 
black hole mass is Mh = 0.95M and its spin parameter is Jh/M"^ = 0.7, with the 
remaining matter forming a torus around the black hole. The subsequent evolution of 
the torus does depend on the strength of the magnetic field. We freeze the spacetime 
metric ( "Cowling approximation" ) and continue to follow the evolution of the torus after 
the black hole has relaxed to quasi-stationary equilibrium. In the absence of magnetic 
fields, the torus settles down following ejection of a small amount of matter due to shock 
heating. When magnetic fields are present, the field lines gradually coUimate along the 
hole's rotation axis. MHD shocks and the MRI generate MHD turbulence in the the torus 
and stochastic accretion onto the central black hole (see Fig.IS]). When the magnetic field 
is strong, a wind is generated in the torus, and the torus undergoes radial oscillations 
that drive episodic accretion onto the hole. These oscillations produce long-wavelength 
gravitational waves potentially detectable by LISA. The final state of magnetorotational 
collapse always consists of a central black hole surrounded by a coUimated magnetic field 
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Figure 5. Magnetoroational collapse of a massive star to a black hole. Snapshots of meridional 
rest-mass density contours, velocity vectors and poloidal magnetic field lines for the initial and 
endpoint configurations for n = 3 collapse. Field lines coincide with contours of vector potential 
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frame denotes the apparent horizon and the shaded region the excision domain. [Adapted from 
|Liu, Shapiro fc Stephens 2007| ] 



and a hot, thick accretion torus. This system is a viable candidate for the central engine 
of a long-soft gamma-ray burst. 



6. Cosmological Growth of Supermassive Black Holes 

Growing evidence indicates that supermassive black holes (SMBHs) with masses in 
the range 10^ — 10^'^ Mq exist and are the engines that power AGNs and quasars. There 
is also ample evidence that SMBHs reside at the centers of many, and perhaps most, 
galaxies, including the Milky Way. The highest redshift of a quasar discovered to date 
is Zqso = 6.43, corresponding to QSO SDSS 1148H-5251 (jFan et al. 2003p . Accordingly, 
if they are the energy sources in quasars (QSOs), the first SMBHs must have formed 
prior to Zqgo — 6.43, or within t = 0.87 Gyr after the Big Bang in the concordance 
ACDM cosmological model. This requirement sets a significant constraint on black hole 
seed formation and growth mechanisms in the early universe. Once formed, black holes 
grow by a combination of mergers and gas accretion. 

The more massive the initial seed, the less time is required for it to grow to SMBH 
scale and the easier it is to have a SMBH in place by Z ^ 6.43. One possible progenitor 
that readilv produces a SMBH is a supermassive star (SMS) with M » lO^M^ ( Rees 
1984 ^hapiro 2004p . SMSs can form when gaseous structures build up sufficient radiation 
pressure to inhibit fragmentation and prevent normal star formation; plausible cosmolog- 
ical scenario s have been proposed that can lead to this situation (jGnedin 2001[ Bromm 
& Loeb 2003) . Alternatively, the seed black holes that later grow to become SMBHs may 
originate from the collapse of Pop III stars < lO'^M© (jMadau fc Rees 200ip . To achieve 
the required growth to ~ IO^Mq by Zq^o ^ 6.43, it may be necessary for gas accretion, 
if restricted by the Eddington limiting luminosity, to occur at relatively low efficiency of 
rest-mass to radiation energy conversion (< 0.2; [Shapiro 2005| and references therein), as 
we discuss below. 

The efficiency of black hole accretion, and the resulting rate of black hole growth, is 
significantly affected by the spin of the black hole. The spin evolution of a black hole 
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m? Characterization 

standard thin disk; nonspinning BH 

turbulent MHD disk 

standard thin disk; photon recapture 

standard thin disk; max spin BH 

Table 1. Rest-mass-to- radiation conversion efficiency vs. black hole spin. 



begins at birth. If the hole arises from the collapse of a massive or supermassive star, 
then it is likely to be born with a spin parmeter in the range ^ a/M < 0.8. The 
lower limit applies if the progenitor star (or core) is nonrotating, the upper limit if it is 
spinning uniformly at the mass-shedding limit at the onset of collapse, as found in the 
simulations of Shibata & Shapiro and Liu & Shapiro 2007 and discussed in Section [5] 
Major mergers with other black holes of comparable mass will cause the black hole to 
spin up suddenly to a/m k, 0.8 — 0.9, as the merged remnant acquires almost all of 
the mass and angular momentum that characterizes a circular orbit, quasistationary 
BHBH binary at the innermost stable circular orbit (ISCO). Once secular gravitational 
radiation loss drives the binary past the ISCO, the black holes undergo a rapid dynamical 
plunge and coalescence, with little additional loss of energy and angular momentum.. This 
anticipated behavior has now been confirmed by numerical simulations of BHBH mergers 
(see Centrella 2007, This Volume, and references therein). These simulations also show 
that gravitational radiation reaction can induce a large kick velocity (> 1000 km/s) in 
the remnants following mergers. While in principle these large kick velocities pose a great 
hazard for the growth of black hole seeds to SMBHs by .Z ~ 6, large kicks are possible 
only if the spins of the black hole binary companions are appreciable and their masses 
are comparable. In the end, gravitational recoil does not pose a significant threat to 
the formation of the SMBH population observed locally, although high mass seeds are 
favored (|Volonteri 2007p . 

Minor mergers with many smaller black holes, isotropically distributed, cause the black 
hole to spi n down: a/m ^ M'"^/^ ([Hughes fc Blandford 20031 Gammie, Shapiro & McK- 
inney 2004). The reason is that the ISCO and specific angular momentum of black holes 
orbiting counter-clockwise is larger than for holes orbiting clockwise, hence the net effect 
of isotropic capture is spindown. 

However, it is likely that most of the mass of a supermassive black hole has been 
acquired by gas accretion, not mergers. Such a conclusion can be inferred from the 
observation that the luminosity density of quasars is roughly 0.1 — 0.2 of the local SMBH 
mass density (jSoltan 1982| lYu fc Tremaine 2002^ , an equality that arises naturally from 
growth via gaseous disk accretion. Steady gas accretion will quickly drive the black hole 
to spin equilibrium, with a spin paranreter that depends on the nature of the flow. If 
accretion occurs via a relativistic "standard thin disk" , then the hole will be driven to the 
Kerr limiting value, a/M = 1 (jBardeen 1970[) . Correcting for the recapture of some of the 
emit ted photons from such a disk reduces the equilibrium spin value to 0.998 ( Thorne 
1974). If, however, the accretion is driven by MRI turbulence in a relativistic MHD 
disk, then recent simulatio ns indicate ([McKinney fc Gammie 2004[ Gammie, Shapiro & 
McKinney 2004' De Vilhc rs et al. 2005^ that the equilibrium spin will fall to ~ 0.95. 
The small differences between these equilibrium spin parameters is deceptive, for they 
correspond to very different rest-mass-to-radiation conversion efficiencies, as shown in 
Table m 
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The significance of these different vaiues is that the growth of a black hole with time 
via steady accretion depends exponentially, on the rest-mass-to-radiation conversion effi- 
ciency, Em, as we will now recall. Define the rest- mass-to-radiation conversion efficiency 
Cm and the luminosity efficiency e^, according to 

eM = L/Moc" = eMia/M) e^ = L/Le (6.1) 

where M is the black hole mass, Mq is the accreted rest-mass, and Le is the Eddington 
luminosity given by 

Le = ^"^^-'"^^ « 1.3 X lO^VeMs erg s'^ . (6.2) 

(tt 

Here r is the growth timescale, 

A/c2 



Le 



OAb^^-^ Gyr , (6.3) 



and /ig is the mean molecular weight per electron. With the above definitions the growth 
rate of a black hole due to accretion is 



dM , , dAfo 



eL(l -Cm) 



fA/ 



T 



Integrating equation (|6.4p for steady accretion with constant efficiencies trivially yields 
the mass amplification of a black hole with time, 



M[t)lM{ti) = exp 



eL(l-eA/)(t-t«) 



Em 



(6.5) 



showing the exponential dependence on the efficiency factors. 

Now a possible clue to the upper limit of e^ is provided by the broad-line quasars 
in a Sloan Digital Sky Survey sample of 12,698 nearby quasars in the redshift interval 
0.1 ^ ^ ^ 2.1. This survey supports the value e^ r^ 1 as a physical upper limit ( McLure 
& Dunlop 2004). Hence to maximize mass amplification via steady accretion, one must 
accrete with a small value of ej\f at the Eddington limit, e^ « 1. 

Fig. [6] evaluates equation (16. 5p for black hole seeds that form at redshift Zi from the 
collapse of Pop III stars in the range 100 < MJMq < 600, and subsequently grow 
to 1O^M0 by Zqso = 6.43. The ACDM concordance cosmological model is assumed. Of 
particular relevance is the case of "merger-assisted" mass amplification, whereby mergers 
account for a tvpic al growth of ~ lO'^ in black hole mass, the remainder being by gas 



accretion (see, e.g.. lYoo fc Miralda-Escude 2004|) . The figure shows that for Zi > 40, the 



required growth of the seed to SMBH status is easily achieved for a rclativistic MHD 
accretion disk, but is only marginally possible for a standard thin disk that accounts for 
photon recapture, and not at all possible for a standard thin disk that drives the black 
hole to maximal spin. However, if the initial black hole seed is less than 6OOM0, accretion 
via a standard thin disk appears to be ruled out altogether. The fact that a black hole 
driven to spin equilibrium by a turbulent MHD disk accretes with low enough efficiency 
to grow to supermassive size by .^qso ~ 6.43 is potentially significant. It points to the 
need for further relativistic MHD simulations of black hole accretion with ever greater 
physical sophistication, including the full effects of radiation transfer. 

These same conclusions also hold if the black hole seed forms much earlier than Zi « 40, 
and they may be tightened if the seed forms later. In fact, it may be likely that the 
seed forms later, at Zi < 40, given that even 4 — cr peaks in the density perturbation 
spectrum for the progenitor halo of SDSS 1148-1-5251 do not collapse until Z ~ 30 in the 
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Figure 6. Black hole accretion mass amplification Mf /Mi versus redshift Zi of the initial seed. 
Here we plot the amplification achieved by redshift Zf = 6.43, the highest known quasar redshift, 
corresponding to 1148+5251. Each solid curve is labeled by the adopted constant radiation 
efficiency, em; the luminosity is assumed to be the Eddington value (ez, = 1). The horizontal 
dashed lines bracket the range of amplification required for accretion alone to grow a seed black 
hole of mass 100 < M/Mq < 600 formed from the collapse of a Pop III star to lO^M©. The 
horizontal dotted lines bracket the required accretion amplification range assuming that mergers 
account for a growth of 10* in black hole mass, the remainder being by gas accretion. [From 
[Shapiro 20051 ] 



ACDM concordance cosmology (see, e.g. Figure 5 in lBarkana fc Loeb 20011 ) Moreover, 
the potential wells of the earliest halos are quite shallow (~ 1 km/s) and may not be able 
to retain enough gas to form stars. Nevertheless, the effect of altering the date of birth 
of the black hole seed is not very great unless Z^ < 20 — 25, as is evident from Fig. 3. 

Should a quasar be discoverd at Zqso substantially above 6.43, it would not be un- 
derstood easily in the context of supermassive black hole growth by gas accretion from 
a seed arising from the collapse of a stellar-mass, Pop III progenitor. 

The analysis presented here is illustrative only; the results may change as the treatment 
is refined. The main point of this example is to emphasize that our understanding of 
structure formation in the early universe as it pertains to the formation and growth 
of supermassive black holes (global physics) depends in part on resolving some of the 
important details of relativistic BHBH recoil and relativistic black hole accretion (local 
physics). To understand these details, in turn, requires large-scale simulations in full 
general relativity, which are now possible and underway. 
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7. Conclusions 

We have discussed a number of black hole scenarios that require numerical simula- 
tions in full general relativity for true understanding. Some of these processes involve 
the formation of black holes, others concern their subsequent growth and interactions. 
The important point is that numerical relativity is at last mature enough to probe these 
issues reliably. Numerical relativity can now serve as an important tool to simulate stellar 
collapse to black holes, the inspiral, merger and recoil of BHBHs, NSNSs and BHNSs, 
the generation of gravitational waves from stellar collapse and binary mergers, accretion 
onto black holes, and countless other phenomena involving black holes and their strong 
gravitational fields. With such a tool, these processes now can be investigated at a fun- 
damental level without many of the ad hoc assumptions and approximations required in 
previous treatments. This should lead to an improved qualitative picture, as well as more 
reliable quantitative results. 

We are grateful to T. Baumgarte, C. Gammie, Y-T Liu, and M. Shibata for useful 
discussions. S.L.S is supported by NSF grants PHY-0205155, PHY-0345151, and PHY- 
0650377 and NASA Grants NNG04GK54G and NNX07AG96G to the University of Illi- 
nois. 
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